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Extreme Wave Events in Directional, Random Oceanic Sea States 



M. Onorato, A. R. Osborne and M. Serio 
Dipartimento di Fisica Generate, Universitd di Torino, Torino, Via P. Giuria, 1 - 10125, ITALY 

We discuss the effects of the directional spreading on the occurrence of extreme wave events. 
We numerically integrate the envelope equation recently proposed by Trulsen et al., Phys of Fluids 
2000, as a weakly nonlinear model for realistic oceanic gravity waves. Initial conditions for numer- 
ical simulations are characterized by the spatial JONSWAP power spectrum for several values of 
the significant wave height, steepness and directional spreading. We show that by increasing the 
directionality of the initial spectrum the appearance of extreme events is notably reduced. 
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We address the occnrrence of extreme wave events 
in numerical models of random, directional, oceanic sea 
states. Extreme ocean waves of this type are of unusually 
large size with respect to the background, surrounding 
waves. They are often refered to as "freak" or "rogue" 
waves and have been given a rough definition in terms 
of an arbitrary threshold, -ffmax > 2.2Hs, for Hg the 
significant wave height 0. Because extreme waves are 
rare, and hence are not often measured , the physical 
mechanisms for their occurrence have not been unequiv- 
ocally established. Possible causes of freak waves in the 
open sea have been investigated by many researchers ||] 
and basically three mechanisms have been proposed: the 
linear interaction of the waves with the currents (geo- 
metrical optic theory) , see ; the simple linear super- 
position of Fourier phases (the resulting large waves are 
also known as transient waves ^) and the modulation 
instability 

The last mechanism and indeed the method addressed 
herein is based on the nonlinear phenomenon known as 
the Benjamin- Feir instability ||^, which describes how 
a monochromatic wave train can be unstable to small, 
side-band perturbations. This physical phenomena has 
been widely studied in wave tank facilities (see for exam- 
ple ||], [0 and references therein) and from a theoret- 
ical and numerical point of view for the 1+1 Nonlinear 
Schroedinger equation (NLS equation in one space and 
one time dimensions ) IlllJl^ and for the fully nonlin- 
ear Euler equations [[13[. A major complication arises 
from the fact the envelope equations (for example NLS) 
are derived from the primitive equations of motion under 
the hypothesis of a narrow-banded spectrum. Higher or- 
der equations in the envelope hierarchy have then been 
proposed ||l^,|l^ in order to overcome the limitation of 
narrow bandwidth. Numerical simulations of these equa- 
tions in 1+1 dimension (l^Jl^ have shown that the prob- 
ability of occurrence of freak waves is strictly related to 
the "enhancement" coefficient 7 and the Phillips param- 
eter a of the JONSWAP power spectrum In [||jl7| 
it was found that increasing the values of 7 and a has 
the effect of increasing the probability of finding freak 
wave events. In addition to the physical limitations as- 
sociated with the narrow spectral width assumption, an 
even more severe limitation of the above results arises 
because these conclusions have been obtained from es- 



sentially one-dimensional numerical simulations. 

In this communication we address the less restricted 
problem of the appearance of freak waves in 2+1 dimen- 
sions, i.e. the wave field depends on the spatial coor- 
dinates X, y and time, t, so that directional spreading 
is systematically included. Preliminary results have al- 
ready been considered in | |T9[ | where simple initial con- 
ditions characterized by a carrier wave plus longitudinal 
and transverse perturbations showed an ubiquitous for- 
mation of large amplitude " freak" waves under the evolu- 
tion of the NLS equation in 2+1. Our aim here, instead, 
is to study the more realistic case characterized by the 
JONSWAP spectrum for a directional wave train, while 
at the same time extending the order of the simulations 
so that energy leakage to high wave numbers does not 
occur (see g). Recently Trulsen et al. proposed a 
model equation (see eqs. (^||) below) that shares many 
of the properties of the Zakharov equation |20j . The main 
advantage of the model is that it is numerically inexpen- 
sive and therefore suitable for extensive numerical simu- 
lations in which large amounts of statistical information 
is desired. The equation can be thought of as a special 
case of the Zakharov equation: 

'^^W^ = y^"Tk,k„k„k3i?(ki,i)*S(k2,t)i3(k3,<) X 

dik + ki - k2 - k3)e'('^+"i--^-'^^'*dkidk2dk3, (1) 

where Tk,ki,k2,k3 = ^(k, ki, k2, ks) is the coupling co- 
efficient (for the analytical form of it see for example 
p2t). We briefiy outline the connection of the Zakharov 
equation with the proposed model. One performs the fol- 
lowing change of variables: k = ko + X with x = (-^i m) 
and 



A{x,t) ^ B{k,t)e- 



-i(w(k)-w(ko))t 



(2) 



For simplicity we perform a simple rotation so that 
kg has the coordinates (A:o,0). The next step con- 
sists in a Taylor's expansion of the coupling coefficient, 
T(k, ki, k2, ka), up to first order in spectral width around 
the carrier wave number, kg. Note that the expansion is 
performed on the nonlinear term, while leaving the linear 
part of the equation unchanged. As a result, the follow- 
ing particular kernel is used in eq. dll) (see ): 
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2fco|Xi-X2| 2fco|Xi-X3| 



(3) 



Under these conditions the triple integral in eq. can 
be performed (see for details) and the following set of 
equations is obtained in physical (configuration) space: 



At+ujaL{d^,dy)A + i 
ojoko 



AV A- 



■-A'Al + ikoA^,l^^ = 



where L — 



,1/4 



1 



.A\'A,- 
(4) 

use of this 



latter operator insures retention of linear dispersion to 
all orders, (j) is the inducead mean flow which is coupled 
with the envelope A via the following nonlocal relation: 



-^F-'{\k\F{\Af}}, 



(5) 



with F the Fourier operator. We call equations (^^ the 
Generalized Dysthe Equation (GDE) since the nonlinear 
part is exactly that of the Dysthe Equation ||l5| , but fur- 
ther includes the linear dispersion relation of the inviscid 
primitive equations of motion. As stated in pl| ], this 
equation shares many properties of the Zakharov equa- 
tion: the four wave resonant interaction curve coincides 
with the Phillips curve and moreover the instability re- 
gion is bounded and is quantitatively similar to the re- 
sults obtained by McLean for the Zakharov equation 
(see Fig. 1 in 1^). 

Since our main objective herein relates to the genera- 
tion of freak waves by the Benjamin-Feir instability and 
since eqs. (|4|-||) show good agreement with the fully non- 
linear inviscid equations for the two dimensional instabil- 
ity diagram, we feel that many of the features of the prim- 
itive equations are retained in the GDE equation. This 
suggests that a study of the generation of freak waves 
using this higher order formulation may indeed allow us 
to investigate many of the effects of directional spreading 
on the generation of freak waves. 

The GDE equations have been solved numerically us- 
ing a standard pseudo-spectral code fully dealiased and 
with second order time-stepping on a doubly-periodic do- 
main with a spatial discretization of 256 x 64. The initial 
conditions consist of a two-dimensional wave field char- 
acterized by the JONSWAP power spectrum. In order to 
statistically search for extreme wave events we focus on 
the temporal evolution of the fourth-order, central mo- 
ment (kurtosis) of the probability density function of the 
wave amplitude. The kurtosis is a statistical measure of 
the probability density function of the wave amplitudes: 
the larger/wider the tails are, the larger are the value of 
the kurtosis. It is clear that when the kurtosis is much 
larger than 3 (the well-known value for a Gaussian pro- 
cess) we are faced, in physical space, with an extreme 



wave event that likely corresponds to a "freak" wave (see 
below for details). 

In order to prepare the initial wave field it is necessary 
to generate a directional, frequency spectrum S{f, 9) — 
P{f)G{e), where P{f) is the original JONSWAP spec- 
trum, and then to transform it into the associated wave 
number spectrum, S{kx,ky). The directional spreading 
function G{9) used here is a cosine-squared function in 
which only the first lobe (relative to the dominant wave 
direction) is considered: 



G{e) 



cos 




2/3^ 



else 



(6) 



(3 is an in important parameter that provides a measure 
of the directional spreading, i.e. as /3 ^ the waves 
become increasingly unidirectional. Using the linear dis- 
persion relation for infinite water depth, uo — g\k.\, the 
wave number spectrum can be written: 



S[kx, ky) — 



|k|3 2|k 



1 -ilhiV °''P 
e Hikij 



(v^-^/fc^)^ 



G{0) 
(7) 

where 6 — arctan(fcj,/fc2,). From eq. (^ the two- 
dimensional surface elevation is computed in the follow- 
ing simple way (see [p5[): 



N M 



■qix, y) = X! X! cos(fcia; + k^ 



(8) 



i=i 3=1 



where the (pij are uniformly distributed random phases 
on the interval (0, 27r) and Cij = yJ^S{ki, kj)AkiAkj. 
Most of our numerical simulations have been started with 
a field of size 5026.5 m x 1885 m, with a dominant wave 
number fco = 0.02 m^^, corresponding to a characteris- 
tic wavelength of ^ 300 m. The time step adopted is 1 
second and the full duration of the simulation for each 
run is one hour in real time. We recall that the nonlin- 
ear time scale Tnl for the effects of the Benjamin-Feir 
instability to appear is of the order of O{l/{e'^ujo)). In 
our numerical simulations the wave steepness is e ^ 0.1 
and the dominant frequency fo can be estimated from 
ko using the linear dispersion relation in infinite water 
depth, fo ~ 0.07 s~^. A rough estimation of the nonlin- 
ear time scale is Tnl — 250 s. For each time step the 
two-dimensional surface has been visualized and the kur- 
tosis has been computed from the space series associated 
with each instant in time of the complete evolution. In 
Fig. [|a, b we show that, for small values of the enhance- 
ment parameter 7 (for 7=1 the significant wave height 
Hg, computed as 4 times the standard deviation of the 
wave amplitude, was estimated to be = 7.26 m ), the 
kurtosis is found to be ~ 3 for all the evolution, indepen- 
dently of the value of the spreading parameter /?. Ad- 
ditional numerical simulations, with many different sets 
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of random phases, have been performed and the results 
are quantitatively similar to the one just presented. A 
more interesting case, from the point of view of extreme 
waves, corresponds to the numerical simulations with ini- 
tial conditions characterized by 7 = 5 {Hg — 9.8 m). The 
main results of these latter simulations arc shown in Fig. 
0a, b. 
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FIG. 1. Kurtosis as a function of time from numerical 
ulations with a = 0.0081 and 7 = 1, for /3=0,10 (a) 
/3=20,30 (b) . 
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ga) indicating the onset of the Benjamin-Fier instability. 
The kurtosis is seen to grow, to reach a maximum value 
and then to decrease and oscillate around the value ~ 3. 
An example of the wave field with the kurtosis — 4.3, 
resulting from the evolution of an initial condition with 
7 = 5 and /? = 10 is given in Fig. |[ a large ampli- 
tude wave with respect to the rest of the quasi-gaussian 
background appears in the middle of the two dimensional 
surface plot; note also that before the large wave a deep 
hole is also evident (this phenomenology is consistent 
with many descriptions of mariners who happened to run 
into a freak wave). As the spreading parameter, /3, is in- 
creased in subsequent numerical simulations (Fig. ^b) 
the kurtosis continues to grow but does not reach values 
as high as in the previous cases. For the larger values of 
/3, the kurtosis oscillates around ~ 3, indicating that the 
probability density function for the wave amplitudes is 
approximately Gaussian. 
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FIG. 2. Kurtosis as a function of time from numerical 
lations with a = 0.015 and 7 = 5, for (3—0,10 (a) and /3 
(b) . 
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Here we see that for small values of the spreading pa- 
rameter, /3, the kurtosis suddenly starts to increase (Fig. 
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FIG. 3. A two dimensional free surface showing a quasi 
gaussian backgraound and a large wave amplitude inthe mid- 
dle of it. The kurtosis of the field is 4.3. The surface 
has been obtained from the nonlinear evolution of a random 
phase wave field characterized by a JONSWAP spectrum with 
a = 0.0081, 7 = 5 and /? = 10. 



A considerable number of numerical simulations have 
been carried out using different values for the selected 
sets of random phases. The results are in qualitative 
agreement with the case just discussed. However, we find 
that the times at which the rare events appear are differ- 
ent for different sets of random phases and that the max- 
imum values of the kurtosis can be also be quite different. 
For a number of random realizations we have found, for 
the case of a quasi-unidirectional spectrum, the appear- 
ance of as many as two extreme events in a single hour of 
simulation time. We note that for very small directional 
spreading (essentially the one-dimensional case) the pat- 
tern of freak wave generation is almost periodic, pri- 
marily because of the well-known phenomenon of Fermi- 
Pasta-Ulam recurrence. In the fully 2-1-1 case the physics 
is much more complicated: as noted in [ p6| , there is a 
permanent frequency downshift which excludes the pos- 
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sibility of recurrence. As a final general consideration it 
should be mentioned that the classical statistical descrip- 
tion of deep water wave trains (using the kinetic equa- 
tions of the Hasselmann, or Boltzmann, type in the quasi- 
Gaussian approximation) is based on the assumption that 
the wave field is spatially homogeneous (see for example 
psf): under this condition the Benjamin- Feir instabil- 
ity is suppressed and the transfer of energy in the power 
spectrum is ruled only by the four-wave resonant inter- 
action whose nonlinear time scale is T/vl = 0(l/(£^/o)). 
Therefore, if the statiscal prediction of the occurence of 
freak waves is desired, a new form of kinetic equation 
should be adopted. This very interesting application is, 
however, beyond the scope of the present paper and will 
be considered in a future work. 
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